/*
 * Triangular.cpp
 *
 *  Created on: Jun 8, 2011
 *      Author: gilang
 */

#include "Triangular.h"

Triangular::Triangular() {
	init();
}

Triangular::Triangular(const int & label) {
	init();
	this->label = label;
}

Triangular::Triangular(Element *e) {
	this->label = e->getLabel();
	this->node = e->getNode();
	this->stiffMatrix = e->getStiffMatrix();
	this->convMatrix = e->getConvMatrix();
	this->convVector = e->getConvVector();
	this->fluxVector = e->getFluxVector();
	this->genVector = e->getGenVector();
	this->material = e->getMaterial();
	this->loadVector = e->getLoadVector();
	this->last = e->isLast();
}

Node & Triangular::operator ()(int i) {
	return node[i];
}

Vector Triangular::constructRHS() {
	return NULL;
}

double Triangular::operator [](int i) {
	return loadVector(i);
}

Matrix Triangular::constructLHS() {
	Vector B=this->B();
	Vector C=this->C();
	Matrix *tempB=new Matrix(TRIANGULAR,TRIANGULAR);
	Matrix *tempC=new Matrix(TRIANGULAR,TRIANGULAR);
	for(int i=1;i<=TRIANGULAR;i++){
		for(int j=1;j<=TRIANGULAR;j++){
			tempB->set(i,j,(getMaterial().getCondX()*B(i)*B(j)));
			tempC->set(i,j,(getMaterial().getCondY()*C(i)*C(j)));
		}
	}
	double c=getMaterial().getWidth()/(2*dimLength());
	Matrix *tempH=NULL;
	double hc=0;
	if(last==true){
		tempH=new Matrix(TRIANGULAR,TRIANGULAR);
		tempH->set(2,2,2);tempH->set(2,3,1);tempH->set(3,2,1);tempH->set(3,3,2);
		hc=getMaterial().getConv()*getMaterial().getWidth()*exposedBoundary(2,3)/6;
		convMatrix+=(*tempH)*hc;
	}
	stiffMatrix+=(*tempB)+(*tempC)+convMatrix;
	delete tempB;delete tempC;delete tempH;
	return stiffMatrix*c;
}

void Triangular::addNode(int i, const Node & n) {
	node[i] = n;
}

double Triangular::dimLength() {
	Matrix mdet = Matrix(TRIANGULAR, TRIANGULAR);
	for (int i = 1; i <= TRIANGULAR; i++) {
		for (int j = 1; j <= TRIANGULAR; j++) {
			if (j == 1) {
				mdet.set(i, j, 1);
			} else if (j == 2) {
				mdet.set(i, j, node[i].getX());
			} else if (j == TRIANGULAR) {
				mdet.set(i, j, node[i].getY());
			} else {
				mdet.set(i, j, 0);
			}
		}
	}

	return mdet.determinant();
}

double Triangular::operator ()(int i, int j) {
	return stiffMatrix(i, j);
}

Vector & Triangular::A() {
	Vector *A = new Vector(TRIANGULAR);
	double a1 = node[2].getX() * node[3].getY() - node[3].getX()
			* node[2].getY();
	A->set(1, a1);
	double a2 = node[3].getX() * node[1].getY() - node[1].getX()
			* node[3].getY();
	A->set(2, a2);
	double a3 = node[1].getX() * node[2].getY() - node[2].getX()
			* node[1].getY();
	A->set(3, a3);
	return *A;
}

Vector & Triangular::B() {
	Vector *B = new Vector(TRIANGULAR);
	double b1 = node[2].getY() - node[3].getY();
	B->set(1, b1);
	double b2 = node[3].getY() - node[1].getY();
	B->set(2, b2);
	double b3 = node[1].getY() - node[2].getY();
	B->set(3, b3);
	return *B;
}

Vector & Triangular::C() {
	Vector *C =new Vector(TRIANGULAR);
	double c1 = node[3].getX() - node[2].getX();
	C->set(1, c1);
	double c2 = node[1].getX() - node[3].getX();
	C->set(2, c2);
	double c3 = node[2].getX() - node[1].getX();
	C->set(3, c3);
	return *C;
}

double Triangular::exposedBoundary(int i,int j)
{
	double lx=node[i].getX()-node[j].getX();
	double ly=node[i].getY()-node[j].getY();
	return sqrt(lx*lx+ly*ly);
}

Triangular::~Triangular() {
	// TODO Auto-generated destructor stub
}
